function res=theta(t,X)
    global b;
    global I;
    global K;
    %res = zeros(1, length(X));
    res(2)=(-b/I)*X(2)+(K/I)*(thetaR(t)-X(1));
    res(1)= X(2);
end